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Abstract 

A procedure is presented for performing three- 
dimensional stress analysis of spiral bevel gears in 
mesh using the finite element method. The procedure 
involves generating a finite element model by solving 
equations that identify tooth surface coordinates. 
Coordinate transformations are used to orientate the 
gear and pinion for gear meshing. Contact boundary 
conditions are simulated with gap elements. A 
solution technique for correct orientation of the gap 
elements is given. Example models and results are 
presented. 

Introduction 

Spiral bevel gears are used to transmit power between 
intersecting shafts. One such application is in 
helicopter transmission systems. In this critical 
application, the gears operate at relatively high 
rotational speed and transmit substantial power (i. e. 
1500 HP at 21,000 rpm). 

Prior research has focused on spiral bevel gear 
geometry 1 ** to reduce vibration and kinematic error, 
improve manufacturability and improve inspection. 
Stress analysis is another important area of ongoing 
research. Accurate prediction of contact stresses and 
tooth root/fillet stresses are important to increase 
reliability and reduce weight. 


Much effort has focused on predicting stresses in gears 
with the finite element method. Most of this work has 
involved parallel axis gears with two dimensional models. 
Only a few researchers have investigated finite element 
analysis of spiral bevel gears 7 * 8 . 

Finite element analysis has been done on a single spiral 
bevel tooth using an assumed contact stress distribution 9 . 
The research reported here will utilize the numerical 
solution for spiral bevel surface geometry to study gear 
meshing. Pinion tooth and gear tooth surfaces will be 
developed based on the gear manufacturing kinematics. 
The individual teeth are then rotated in space to create a 
multi-tooth model (4 gear and 3 pinion teeth). The tooth 
pair contact zones are modeled with gap elements. The 
model development procedure and finite element results 
are presented. 

Equations for Tooth Surface Coordinates 

The system of equations, briefly summarized here, re- 
required to define the coordinates of a face-milled spiral 
bevel gear surface were developed by Handschuh and 
Litvin 9 . 

A conical cutting head, attached to a rotating cradle, 
swings through the work piece. Parameters U and 6 
locate a point on the cutting head in coordinate system S c 
attached to the cutting head as shown in figure 1 and 
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described by the following equations. 

T cotij/ - U cost}/ 
U sinty sin0 
c U sinij/ cos0 

1 


( 1 ) 


The roll angle of the cradle <f> c is used to locate the 
rotating cradle with respect to the fixed machine 
coordinate system S m . Parameters U, 9 and <f> c9 along 
with various machine tool settings can be used to 
completely define the location of a point on the 
cutting head in space. 

Since the kinematic motion of cutting a gear is 
equivalent to the cutting head meshing with a 
simulated crown gear, an equation of meshing can be 
written in terms of a point on the cutting head (i.e., 
in terms of U, 9 and The equation of meshing 
for straight-sided cutters with a constant ratio of roll 
between the cutter and work piece is given in 1,5,6 as: 

(U - r cot yf/ cos yf/) cos y sin r + S(m cw - sin 7) 
cos \f/ sin 9 T cos 7 sin 7 sin (q - 
± E m (cos 7 sin yf/ + sin 7 cos \p cos r) 

- L* sin 7 cos yf/ sin r = 0 (2) 

The upper and lower signs are for left and right hand 
gears respectively. The following machine tool 
settings are defined 5,6,9 . 

yf/ cutting blade angle 

r (H q ± 

q cradle angle 

7 root angle of work piece 

E m machining offset 

vector sum of change of machine 
center to back and the sliding base 
m^ 4>J4 > w , the relationship between 

the cradle and work piece for a 
constant ratio of roll 

U generating cone surface coordinate 

S radial location of cutting head in 

coordinate system S m 
r radius of generating cone surface 


coordinates of a spiral bevel gear. The 3 parameters U, 
9 and <f> c are defined relative to the cutting head and cradle 
coordinate systems (S c and SJ respectively. These 
parameters can be transformed through a series of 
coordinate transformations to a coordinate system attached 
to the work piece. Or U, 9 , <f> e can be mapped into X„, 
Y w , Z* in coordinate system S w attached to the work 
piece. These transformations, used in conjunction with 
two other geometric requirements, give the two additional 
equations. 

The correct U, 9 and <f> e that solves the equation of 
meshing, must also, upon transformation to the work 
piece coordinate system S w , result in a axial coordinate Z* 
that matches with a preselected axial position Z. (See 
figure 2) 


Z w - Z = 0 (4) 

This equation along with the correct coordinate 
transformations (see equation 11) result in a second 
equation of the form 

f 2 (U 9 6, 4> c ) =0 (5) 

A similar requirements for the radial location of a point 
on the work piece results in (See figure 2). 

r -( xl * = 0 ( fi ) 


The appropriate coordinate transformations (see equation 
11) will convert equation (6) into a function of U, 9 , and 

<t> c - 


/ 3 «/,e,< |> e ) = 0 (7) 

Equations (3), (5) and (7) form the system of nonlinear 
equations necessary to define a point on the tooth surface. 


Solution Technique 

An initial guess U°, 9 °, <£ c ° is used to start iterative solu- 
tion procedures. Newton’ s method is used to determine 
subsequent values of the updated vector (TJ k , 0 k , ^ c k )‘ ,0 . 


This is equivalent to: 

/, ( u , 0, = 0 (3) 

U k 

0* 


U 

e*-i 

► + < 

Yf' 1 

Because there are 3 unknowns U, 9, and <f> c ; three 
equations must be developed to solve for the surface 

< 


.*-1 
T( J 


jf 1 , 
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Where the vector Y is the solution of: 


dU 


dU 


du 



a/ 2 (4>!'V 

dU 

d$ c 

a/ 2 (e*-‘) 

a/ 2 (0 

du 

a*. 


a/,<4>r') 

du 

34> t 

fiiu 1 ' 1 . 

e»-‘. *r> 


e*-\ 

AW*- 1 , 



yf - 1 

if 1 


( 9 ) 


The 3x3 matrix in the preceding equation is the 
Jacobian matrix and must be inverted each iteration 
to solve for the Y vector. The equation of meshing, 
function f, is numerically differentiated directly to 
find the terms for the Jacobian matrix. Function f 2 
and f 3 can not be directly differentiated with respect 
to U, 8 and tj> c . After each iteration U k '\ 8 k '\ <t> kA 
(in the cutting head coordinate system Sc) are 
transformed into the work piece coordinate system, 
S w , with the series of coordinate transformations as 
given in equation 1 1 . 



r cot ijr - U cos \|r 

U sin \|/ sin 0 

U sin cos 0 


( 10 ) 


numerical differentiation into ^ + inc, Y w + inc, Z* + 
be. Equations (4) and (6) are then used to evaluate the 
numerical differentiation. Function f,, f 2 , f 3 and the 
partial derivatives of f,, f 2 , f 3 required to the Jacobian 
matrix are updated each iteration. The iteration contbues 
until the Y vector is less than a predetermined tolerance. 
This completes the solution techniques for a sbgle pobt 
on the spiral bevel gear surface. 

The four comers of the active profile are identified from 
the tooth geometry plane as shown b figure 2. Pobt 1 
on the surface is chosen to be the lowest pobt of the 
active profile on the toe end. The initial guess to start the 
procedure has to be sufficiently close to the correct 
solution for convergence to occur. The solution for the 
first pobt proves to be an adequate initial guess for any 
subsequent pobts on the surface. 

Subsequent bterior surface pobts are found by 
berementbg r = (X 2 + Y 2 ) 1 ' 2 and Z. By adjustbg the 
berements used, a surface mesh of any density can be 
calculated. The process is repeated four times for each of 
four surfaces; gear convex, gear concave, pinion convex 
and pinion concave. 

Sbce all four surfaces are generated bdependently, 
additional matrix transformations are required to obtab 
the appropriate orientation for meshbg. The proper 
convex and concave surface orientation (for both the gear 
and pinion) is found by fixbg the concave surface and 
rotatbg the convex surface until the correct tooth 
thickness is obtabed. The correct angle of rotation is 
obtabed by matchbg the tooth top land thickness with the 
desired value. 


where: 


["j = 

[*„/«►.)] wj [M p J [M„/( <U [MJ (11) 


Each matrix [M] above represents a transformation 
from one coordbate system to another. (See 
Appendix I for the specific matrices.) 

Functions f 2 and f 3 are evaluated by startbg with an 
initial U k , 0 k and <f> k , performbg the transformation 
b equation 11 and evaluatbg equations (4) and (6). 
The numerical differentiation of f 2 and f 3 is 
performed by transformbg U k + be, ^ 4- be, <j> k 
+ be (where be is a small berement appropriate for 


Gear and Pinion Orientations Required For Meshing 

After generatbg the pinion and gear surface as described 
above, the pinion cone and gear cone apex will meet at 
the same pobt as shown b figure 3. This pobt is the 
origb of the fixed coordbate system attached to the work 
piece bebg generated. To place the gear and pinion b 
mesh with each other rotations described b the followbg 
example are required: 

1. The 19th tooth is selected for meshbg. This 
corresponds to rotatbg the generated tooth 190° CW 
about the Z* axis. (For this example, N t = 36 total teeth 
on the gear.) For the general case the gear tooth is 
rotated 360/N t +180° CW about its axis of rotation (ZJ. 
This corresponds to selectbg the ith tooth of the gear to 
be b mesh where i = N t /2 + 1. 

2. The pinion is rotated by 90° CCW about the Y axis. 
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Note: rotation (1) corresponds to selecting a different 
tooth on the gear to be in mesh; however, rotation (2) 
corresponds to physically rotating the entire pinion 
until it meshes with the gear. 


This vector is written relative to S c> a coordinate system 
fixed to the cutting head. To obtain the vector n w relative 
to the coordinate system fixed to the work piece, the 
following transformation must be performed. 


3. Because the four surfaces are derived 
independently, their orientation is random with 
respect to meshing. The physical location of the gear 
and pinion after rotations (1) and (2) correspond to 
the gear and pinion in mesh with severe interference. 
To correct the interference the pinion is rotated CW 
about its axis of rotation (ZJ until surface contact 
occurs. For this example the rotation was 3.56°. 
Figure 4 shows an example of a simulated gear pair 
meshing. The generated pinion tooth was copied and 
rotated 12 times and the generated gear tooth was 
copied and rotated 36 times. 

Contact Simulation 

The tooth pair mesh contact point can be located by 
a method described by Litvin 6 or by a search 
technique. Pairs of finite element node points (one 
on each tooth surface) are evaluated until the pair 
with smallest separation distance is obtained. (A 
finer finite element mesh would yield greater 
resolution.) Once the contact point is established, a 
vector normal to the surface at the contact point is 
calculated. 

The intersection of the normal vector on the pinion at 
the contact point with the gear surface identifies the 
second point required to define the gap element. 
Additional gap elements are obtained by taking 
additional vectors from other pinion surface finite 
element nodal points (parallel to the contact point 
normal vector), and calculating where they intersect 
the mating gear surface. Finite element nodal points 
on the gear surface are located to the intersection 
points of the normal vectors and the gear tooth 
surface. 

The vector normal to the cutting surface is given in 5 
as: 


(13) 

Where [L*J is found by removing the 4th row and 
column from [M^J. 

Each gap element is connected between a node on the 
pinion surface, hereafter designated node 1, and the 
corresponding intersection node point on the gear surface, 
designated node 2. The intersection point on the gear 
surface is found as follows: 

Consider node 1, a point on the pinion with coordinates 
P x , Py , P x in S w . Let Q*, Q y , Q z be any point in space 
such that (Q*, Q y , QJ l = bfo, ziy, nj* where b is a scale 
factor, and n*, iiy and n, are the components of the normal 
vector at the contact point. The intersection of the normal 
vector from node 1 with the mating gear surface defines 
node 2 and has to satisfy the following three equations: 

Q* =G, = X,, + b n„ 

Qy = G y = X,, + b iiy (14) 

Qx =0, = Zp + bn, 

Where G x , G y , G z is a point on the gear surface. A point 
on the gear surface must also satisfy 


(G 1 


r cos ^ * e ~ U cos i|r c 

X 

G y 

= [MJ 

Us in i|r c sin0 

G 

U sin i|r c cos 0 



1 


Equation IS leads to a system of 3 nonlinear algebraic 
equations. These three equations, along with the equation 
of meshing for the gear surface provide a system of 4 
equations and 4 unknowns (u, 6 , <f> c , b). these equations 
are solved with Newton’s method described earlier. The 
intersection of the normal from node 1 on the pinion with 
the gear surface is now obtained. 


{"»} 


sin 

cos sin x 

cos cos x 


( 12 ) 


This procedure is used to locate the intersection of 
normals from all points on the pinion surface (in the 
contact zone) with the gear surface. The gear tooth 
surface is remeshed utilizing the intersection points as 
shown in figure 5. Gap elements are connected between 
corresponding nodal points on the pinion and the 
intersection points on the gear surface. 
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Finite Element Model 

An example analysis was performed using gears from 
the NASA Lewis Spiral Bevel Gear Test Facility. In 
this case, the left hand pinion mates with the right 
hand gear. Counter clockwise rotation of the pinion 
results in the concave surface on the pinion mating 
with the convex gear surface. 

The design data for the pinion and gear are given in 
Table I. The design data are used with methods 
given in 6 to determine the machine tool settings for 
the straight sided cutter data given in Table II. 

The finite element gear pair model, shown in figure 
6, contains 4 gear teeth and three pinion teeth. The 
model had 10,101 nodes (30,303 degrees of freedom) 
and 7596 eight noded 3D brick elements. 

The three pinion teeth are fixed in space with zero 
displacements and rotations. This is done by setting 
x, y and z displacements equal to zero on the four 
comer nodes of each rim section. The gear is 
constrained to rotation about its axis of rotation. The 
gear is loaded with a torque of 9450 in lbs on the 
gear by applying 4725 lb force located 2.0976" from 
the gear axis of rotation. 

At the orientation chosen between the pinion and gear 
two pairs of teeth were in contact. One pair had 
contact near the middle region of the tooth and 
another pair had contact near the toe (i.e. , about to 
go out of mesh). Initially the model started with a 
total of twenty one gap elements. For the tooth that 
is approximately midway through mesh, fifteen gap 
elements were used. For the tooth about to leave 
mesh, four gap elements were used. The analysis 
starts with one gap elements closed in each contact 
zone. Within the finite element code an iterative 
process is used to determine how many gap elements 
must close to reach static equilibrium. The solution 
iterated four times before reaching equilibrium with 
four gap element closed in the main contact region 
and one closed in the edge contact region. Stress 
contours for the two pinion teeth with contact are 
shown in Figures 7 and 8. 

The average nodal minimum principal stresses in the 
main contact zone average -204,000 psi with a 
maximum of -299,900 psi and a minimum of 
-103,574 psi. The corresponding elemental stresses 
average -103,500 psi with a maximum and minimum 
of -123,900 psi and -79,500 psi, respectively. 
Thenodal stresses are higher because of load 


concentration from the gap elements. These stress ranges 
bracket the estimated hertzian stresses for the gear set 
under the same load conditions. Contact with only 4 gap 
elements, along with large stress gradient among adjacent 
nodes indicate the need for a finer finite element mesh for 
improved stress prediction. 

Conclusions 

A multi tooth finite element model was used to perform 
three-dimensional contact analysis of spiral bevel gears in 
mesh. Four gear teeth and three pinion teeth are 
generated by solving the equations, based on gear 
manufacturing kinematics, that identify tooth surface 
coordinates. The gear and pinion are orientated for 
meshing with coordinate transformations. 

Surface stresses are evaluated with gap elements. The 
gear surface is remeshed with nodal points relocated to 
identify the correct gap element orientation. Initial FEA 
stress results compare favorably with calculated hertzian 
contact stresses. However, large stress gradients between 
adjacent nodes in the contact zones indicate a need for 
greater finite element mesh refinement. 
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Appendix I 




10 0 0 

0 cos <{> e Tsin 4> e 0 
0 ±sin<J> e cos<t> c 0 

0 0 0 1 


( 17 ) 


Matrix [MJ transforms coordinate system S m into system 
S p which orientates the apex of the gear being 
manufactured. 


HV = 


cos 6 0 -sin 6 -L m sin 6 
0 10 ±E m 

sin 6 0 cos 6 L m cosb 
0 0 0 1 


(18) 


Matrix [MJ transforms coordinate system S p into system 
S # which located the apex of the gear being manufacturer 
with respect to S m 


Coordinate transformations involving both rotation 
and translation require mixed matrix operations of 
multiplication and addition. Matrix representation of 
coordinate transformations will need only 

multiplication of matrices if position vectors are 
represented by homogeneous coordinates'. The 

following4 x 4 matrices are required to transform the 
homogeneous coordinates of a point on the cutting 
head to a point on the work piece. 

Matrix [MJ transforms coordinate system S c , 
attached to the cutting head, into system S„ rigidly 
connected to the cradle. 


WJ « 


cosp 0 sinp 0 

cos<t> w 0 0 

-sinp 0 cos p 0 

0 0 0 1 


(19) 


Matrix [M,J transforms coordinate Sa into system Sw 
attached to the work piece. 


. twj = 


coscj)^ ±sin<(> w 0 0 

tsih«|> k cos4> w 0 0 


( 20 ) 


0 0 10 
0 0 0 1 


IMJ 


10 0 0 
0 cos q Tsin q *sin q 
0 ±smq cosq scosq 
0 0 0 1 


(16) 


Matrix [MJ transforms coordinate system S, into 
system S m attached to the frame. 


TABLE I - PINION AND GEAR DESIGN DATA 



PINION 

GEAR 

Number of teeth pinion 

12 

36 

Dedendum angle, deg 

1.5666 

3.8833 

Addendum angle, deg 

3.8833 

1.5666 

Pitch angle, deg 

18.4333 

71.5666 

Shaft angle, deg 

90.0 

90.0 

Mean spiral angle, deg 

35.0 

35.0 

Face width, mm (in) 
Mean cone 

25.4(1.0) 

25.4 (1.0) 

distance, mm fin) 
Inside radius of gear 

81.05(3.191) 

81.05 (3.191) 

blank, mm (in) 

5.3 (0.6094) 

3.0 (.3449) 

Top land thickness, mm (in) 

2.032 (0.080) 

2.489 (.098) 

Clearance, mm (in) 

0.762 (0.030) 

0.92964 (0.0366) 
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TABLE II - GENERATION MACHINE SETTINGS 


PINION 


GEAR 


Radius of cutter, r, in 
Blade angle, deg 
Vector sum, L* 

Machine offset, E m , in 
Cradle to cutter distance. 
Cradle angle q, deg 
Ratio of roll, M w 
Initial cutter length, u, in 



Concave 

Convex 

Concave 

Convex 


2.9656 

3.0713 

3.0325 

2.9675 


161.9543 

24.33741 

58.0 

22.0 


0.038499 

-0.051814 

0.0 

0.0 


0.154576 

-0.1742616 

0.0 

0.0 

, in 

2.947802 

2.8010495 

2.85995 

2.85995 


63.94 

53.926 

59.234203 

59.234203 


0.30838513 

0.32204285 

0.950864 

0.950864 


9.59703 

7.42534 

8.12602 

7.89156 

, deg 

126.83544 

124.43689 

233.9899 

234.9545 

•« deg 

-0.85813 

-11.38663 

-0.35063 

-12.3384 


*c 



Figure 1. Cutting head cone surfaces and attached coordinate system. 
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Figure 3. Orientation of gear and pinion based on solution of the system of equations, and after 
rotations required for mesh. 
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7 X 

7 

Figure 4. Complete 3D model of gear and pinion in mesh. 



Figure 5. Distorted gear after connecting gap elements. 
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Figure 6. Seven tooth model used for finite element analysis of mesh. 
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Figure 7. Stress contours in pinion tooth. 
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